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1.  Introduction 

There  are  several  aspects  of  many  observed  univariate  time  series  which 
are  not  satisfactorily  accounted  for  in  standard  time  series  analysis:  they 
include  non-Gaussian  marginal  distributions,  dependence  beyond  second-order 
moments  (autocorrelations)  and  directionality  in  the  time  series.  Quite 
often  a  Gaussian  distribution  will  be  inappropriate  because  the  variable 
being  modelled  has  a  positive  and  highly  skewed  distribution,  i.e.,  the 
service  times  in  a  queue  or  daily  flows  in  a  river.  Many  particular  such 
distributions  can  be  envisaged  and  time  series  models  have  been  constructed 
for  them.  Examples  are  Gamma  distributions  (Gaver  and  Lewis,  1980;  Lewis, 
1981;  McKenzie,  1982;  Lawrance,  1982)  and  mixed  exponential  distributions 
(Gaver  and  Lewis,  1980;  Lawrance,  1980a;  Lawrance  and  Lewis,  1982) 

However  the  simplest,  most  widely  used  and  most  tractable  analytically 
of  these  distribution  models  is  the  exponential  distribution.  Like  Gaussian 
random  variables,  exponentially  distributed  random  variables  enjoy  many 
special  properties;  also  they  can  be  mildly  transformed  quite  easily  into 
distributions  which  are  more  skewed,  or  less  skewed, than  the  exponential.  The 
Weibull  distribution  is  an  example,  being  just  a  power  transformation  of  an 
exponentially  distributed  random  variable.  Thus  the  approach  here,  following 
earlier  work  (Gaver  and  Lewis,  1980;  Lawrance  and  Lewis,  1980,  1981)  is  to 
regard  the  exponential  variables  as  canonical  and  to  develop  their  use  in 
time  series  modelling.  It  should  also  be  noted  that  time  series  of  uniformly 
distributed  random  variables  are  obtained  by  exponential  transformations  of 
the  exponential  time  series;  this  has  many  possibilities  and  the  process 
of  uniformly  distributed  random  variables  will  be  considered  elsewhere. 

The  work  cited  above  has  concentrated  for  the  most  part  on  first-order, 
non-Gaussian  autoregressive  models,  both  of  the  standard  type  (constant 
coefficient  additive  linear  combinations)  and  a  random  coefficient  type 
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introduced  by  the  authors.  The  extension  of  the  models  to  higher  order 
autoregression  is  clearly  necessary  to  attain  flexibility  in  modelling  cor- 
relation and  dependency  structure  of  the  processes,  but  these  extensions  are 
in  no  way  as  immediate  as  in  the  standard  linear  Gaussian  case.  A  simple 
mixing  device  can  be  used  (Jacobs  and  Lewis,  1983)  but  the  range  of  correla- 
tions attained  is  much  narrower  than  the  range  attained  in  the  standard 
linear,  second-order  autoregressive  structure.  A  broader  extension,  called 
the  EAR(2)  model,  was  obtained  in  the  exponential  case  by  Lawrance  and  Lewis 
(1980),  but  its  innovation  variable  has  a  zero  component  and  this  will  often 
be  hard  to  justify. 

A  major  part  of  the  present  work  consists  of  obtaining  a  very   broad  and 
rich  extension  of  the  NEAR(l)  model  (Lawrance  and  Lewis,  1981)  to  a  second- 
order  autoregressive  process;  this  NEAR(2)  model  was  first  proposed  in 
Lawrance  (1980b)  and  later  reviewed  in  Raftery  (1981),  although  necessary 
analysis  of  its  innovation  structure  was  not  given.  The  innovation  random 
variable  for  the  NEAR(2)  process  is  proved  here  to  exist  without  unnatural 
boundaries  on  its  (four)  parameter  region  and  does  not  have  a  zero  component. 
Taken  out  of  the  context  of  second-order  autoregressive  exponential  time 
series  this  result  gives  a  very   broad  structure  for  combination  of  exponential 
random  variables.  In  fact  it  gives  a  way  to  combine  three  (possibly  dependent) 
exponential  random  variables  with  an  independent  pair  of  (possibly  dependent) 
exponential  random  variables  to  give  another  exponential  random  variable.  The 
result  has  a  number  of  uses  beyond  the  second-order  time  series  model  giving, 
for  instance,  a  first-order  model  which  is  broader  than  the  NEAR(l)  model,  and 
several  second-order  mixed  autoregressive-moving  average  models.  Schemes  for 
obtaining  negative  correlation  can  also  be  accommodated.  These  models,  as  in 
the  first-order  case,  include  aspects  of  dependence  beyond  autocorrelation  and 
also  directionality  in  the  time  series. 


The  richness  of  the  four-parameter  NEAR(2)  model  and  the  fact  that  an 

infinite  number  of  cases  of  the  model  with  identical  correlation  structure 

are  available  forces  consideration  of  higher-order  aspects  of  dependence. 

The  analysis  of  the  higher  order  aspects  of  the  exponential  time  series  is 

at  a  fairly  early  stage  and  is  as  follows.  First  it  will  be  shown  that  the 

autocorrelations  p{i),   z   =  0,  +  1,  +2,...  for  the  NEAR(2)  process  satisfy 

the  Yule-Walker  equations  with  constants  a-,  and  a?  which  are  functions 

of  the  four  parameters  of  the  model.  This  follows  immediately  from  the  fact 

that  X       is  a  random-coefficient,  linear  additive  combination  of  X  ,  , 
n  n-1  ' 

Xn-2  and  the  innovatlon  r-v-  En  •  Secondly,  it  can  be  shown  (Lewis  and 
Lawrance,  1983a)  that  the  residuals  X  -  a,X  ,  -  a?X  2  ,  which  are  the 
usual  residuals  for  second-order  constant  coefficient,  linear  additive 
autoregressive  processes  are  uncorrelated. 

Thus  although  the  standard  analysis  of  time  series  stops  with  uncorrelated 
residuals,  i.e.  a  flat  spectrum  for  the  residuals,  such  residuals  can  also  be 
used  to  good  effect  to  investigate  higher  order  aspects  of  dependence  in  the 
NEAR(2)  model.  In  fact,  if  the  autoregression  is  not  of  the  standard  type 
(constant  coefficient,  linear  additive)  the  (uncorrelated)  residuals  will 
have  non-random  scatters  and  higher  order  dependence  properties.  In  partic- 
ular, the  squared  residuals  will  have  non-zero  autocorrelations  and  the 
cross-correlations  of  residuals  and  squared  residuals  except  in  the  zero- 
lag  case  will  be  non-zero;  both  sets  of  correlations  are  effectively  zero 
when  a  standard  second-order  autoregressive  model  is  appropriate. 

An  extension  of  this  residual  analysis  is  suggested  by  the  fact  that 

autocorrelations  are  time  reversible,  while  the  process  itself  is  not  time 

reversible,  i.e.  is  directional.  Thus  one  constructs  reversed  residuals 

X  -  a,X  ,,   -  a0X  ,0  which  are  to  be  used  to  assess  directionality.  Analysis 
n    1  n+1    2  n+2  J  J 


of  both  types  of  residuals  can  be  fruitful;  because  of  space  limitations  we 
present  here  only  an  investigation  of  cross-correlations  of  residuals  and 
squared-residuals.  This  particular  case  of  the  possible  residual  analyses 
is  illustrated  by  some  theoretical  calculations,  some  simulations  and  a  brief 
application  to  a  transformed  series  of  wind  speed  data.  The  complete  residual 
analysis  including  directional  aspects  and  their  investigation  using  the 
additional  tool  of  reversed  residuals  will  be  given  elsewhere. 


2.  Exponential  Time  Series  Models 

Our  aim  in  this  section  is  to  give  in  outline  the  ideas  leading  to  the 
time  series  model  of  main  concern  in  this  paper,  and  called  NEAR(2),  follow- 
ing the  earlier  terminology  NEAR(l).  The  NEAR(2)  model  has  four-parameters, 
and  incorporates  and  broadens  the  earlier  two-parameter  EAR(2)  model  (Lawrance 
and  Lewis,  1980).  The  NEAR(2)  model  will  be  exponential  in  marginal  distri- 
bution, have  second-order  autoregressive  Markov  dependence,  and  have  auto- 
correlations satisfying  second-order  difference  equations  of  the  familiar 
Yule-Walker  type.  In  addition  it  will  have  dependence  beyond  autocorrelation, 
and  will  have  adjustable  directionality.  It  is  not  linear  in  the  standard 
sense,  having  random-coefficient,  linear  additive  autoregressive  structure, 
but  neither  is  it  non-linear  in  the  standard  sense  of  incorporating  powers  or 
products  of  lagged  variables. 

Writing  {X  }  for  the  time  series  variables,  and  {E  }  for  an  i.i.d. 
exponential  innovation  sequence  of  unit  mean,  the  two-parameter  NEAR(l)  model, 
as  previously  defined,  is  given  by 


Xn  = 


exn  ,    W.p.  a 


Yi-1 


'  E    w.p.  p  , 


n 


w.p.  1-a 


+  <  (2.1) 

bEn    w.p.  1-p  , 


with  b  =  (l-a)e  and  p  =  (l-$)/{l  -  (l-a)3)  .  The  parameter  region  is  in 
general  0<_a,  B_<l,a=B?M  .  The  case  3  =  1 ,  0  <_  a  <  1  is  rather 
special,  and  has  been  called  the  TEAR(l)  model,  and  when  a  =  1,  0  <_  $  <  1  , 
the  earlier  EAR(l)  model  is  recovered.  Except  for  this  EAR(l)  case,  the 
NEAR(1)  model  does  not  allow  zero  innovations  (Gaver  and  Lewis,  1980)  and  so 
is  more  statistically  acceptable.  The  zero  innovation  implies  that 
X  =  3X    and  thus  3  can  be  determined  exactly  from  runs  down  in  the 
sample  path  of  the  process. 

In  general  the  i.i.d.  innovation  in  the  NEAR(l)  process  are  formed  as 
the  probabilistic  mixture  of  two  exponentials,  and  are  thus  easily  simulated. 


The  NEAR(2)  model  is  a  direct  generalization  of  (2.1)  and  takes  the  form 


VI 


B1Xn_1  W.p.    0] 

B2Xn-2         w,p'   a2 


+  E 


(2.2) 


W.p.  1-ou  -  a? 


with  parameter  region  0  £  s]  ,  e2  <  1 ,  a]  >  0,  a2  >  0,  a]  +  a2  <  1  ;  (E  >  is 
an  appropriately  chosen  innovation  sequence.  Many  special  cases  can  arise 
when  the  above  restrictions  include  some  of  the  equalities  and,  for  the  pur- 
poses of  a  general  development,  it  is  best  to  regard  the  inequalities  as 
strict.  Given  that  {XR}  is  required  to  have  an  exponential  marginal  distri- 
bution, the  main  question  concerns  whether  there  is  a  valid  probability 
distribution  for  En  .  The  Theorem  proved  in  Section  2.3  will  show  that 
this  is  the  case,  and  that  the  distribution,  when  the  inequalities  on 
a.|  f  a2  and  ^  ,  e^  -jn  the  parameter  region  are  strict,  takes  the  form 


E  = 


b2En 


b^E 
3  n 


w.p.  l-p9  -  p„  , 


2   K3 


w.p.  p9  , 


w.p.  p 


3  ' 


(2.3) 


a  probabilistic  mixture  of  three  exponentials.  To  establish  this  result  a 
fairly  detailed  analysis  of  a  derived  moment  generating  function  is  required. 
This  is  necessary  since  a  direct  moment  generating  function  solution  of  (2.2) 
for  E   does  not  establish  that  E   has  a  proper  distribution;  all  that  is 
shown  is  that  the  solution  is  a  possibly  improper  mixture  of  three  exponentials. 


3.  Validity  of  the  NEAR(2)  Model 

In  this  section  we  prove  the  following 
THEOREM.  Let  {E  }  be  an  i.i.d.  sequence  of  unit  mean  exponential  random 
variables.  Then  if  the  four  parameters  a,,  &«»  (?■. ,  B2  satisfy 
a,  >  0,  a2  >  0,  a,  +  a2  <  1 ,  0  <  3, ,  3  <  1 ,  the  relationship 


6lXn-l    w,p'  al  ' 

Xn  =  i32Xn-2    w'p'  a2  ' 

0         W.p.  1-a,  -  a? 


+  En  >   n  =  0,  +  1,  +  2,...,    (3.1) 


where 


En  "  b2En 


biE„ 
3  n 


w.p. 
w.p, 
w.p. 


l-p2  -  P3  , 
P2  » 
p3  ' 


(3.2) 


defines  a  stationary  sequence  of  (marginally)  exponentially  distributed 
random  variables  with  mean  one.  Here 


and 


P2  = 
p3  = 

0  <  b„  = 


{(a]31    +  a232)b2   -    (a-,    +  a2)&}  32>/{  (t>2   "   b3)0    "    t>2)} 
{(c^    +  02)3^2  "    (a131    +   a232)b3)/{(b2  -   b3)(l    -  b3)} 


where 


S   =    (1    -  ■a1)B1    +   (1    -  a2)32    , 


and 


(3.3) 
(3.4) 


1  {s  -   (s2  -  4r)1/2}  <  b2  =  1  |s  +  (s2  -  4r)}  <  1    ,  (3.5) 


(3.6) 


r  =    (1    -   a-j    -   a2)3-|32    • 


(3.7) 


PROOF.  For  the  NEAR(2)  model  specified  by  (3.1)  -  (3.7),  let  <|>Y(t)  and 

A 

4>E(t)  be  the  moment  generating  functions  of  the  {X  }  and  {E   }  sequences; 
then  quite  simply  if  stationarity  of  the  {X^}  series  is  assumed, 

<j>x(t)  -  4>E(t){a1*x(B1t)  +  a2<}.x(62t)  +  (1  -  ci-j  -  a2)}  .  (3.8) 


Assuming  an  exponential  marginal  distribution  of  unit  mean  for  {X  }  ,  then 
the  independent  distribution  of  £E  }  has  moment  generating  function, 
possibly  not  proper,  given  by 


(l+31t)(l+6?t) 

*E(t)  =  -Y " •      (3.9) 

(l+t)[(l-a1  -  a2)6]32t  +  {(1-a^B-,  +  (l-a^K  +  1] 


First  we  remark  that  (3.9)  simplifies  in  special  cases,  such  as  when  one  of 
3^  >  32  =  1 ,  3-j  =  32>  or  a^  +  a2  =  1  .  For  simplicity  in  the  development, 
the  parameter  region  considered  in  the  theorem  is  that  defined  by  0  <  3, , 
$2  <  1 ,  a-j ,  a2  >  0,  a,  +  a2  <   1  .  Interesting  special  cases  are  more  easily 
treated  separately;  these  include  NEAR(l)  models  and  particularly  tractable 
types  of  NEAR(2)  models. 

It  is  convenient  to  establish  right  away  that  the  quadratic  term  in  the 
denominator  of  (3.9)  has  real  distinct  and  positive  roots,  b-,  and  b2  ;  this 
eliminates  any  subsequent  need  to  invert  such  a  term  as  a  whole.  The  required 
condition  for  real  distinct  roots  is  that 


2 

[(l-a-,)31  +  0-a2)B2]  "  4(l-a-,  -  a^g^ 


be  positive:  this  is  so  from  its  equality  to  the  expression 


[(l-a1)3]  -  (l-a2)e2]   +  4a1a23132 


which  is  clearly  positive;  the  positivity  of  the  roots  b,  and  b?  is  obvious 
from  (3.9)  since  their  product  and  sum  given  in  (3.11)  and  (3.12)  below  are 
both  positive. 

Subject  to  the  above  qualifications  concerning  b-,  and  b?  ,  a  partial 
fraction  expansion  of  (3.2)  can  be  written  in  the  suggestive  form 

♦  E(t)  -  <1-p2-p3)  -^  ♦  p2  j^  *  p3  t^-z  .  (3.10) 

Comparisons  between  (3.9)  and  (3.10)  then  show  that  b2,  b3  and  p2,  p3  may 
be  obtained  in  terms  of  3-,,  $2  and  a-,,  a2  by  solving  the  equations 

b2  +  b3  =  (l-a])31  +  (l-a2)62  ;              (3.11) 

b2b3  =  (1-a-,  -  a2)$132  ;  (3.12) 

(l-b2)p2  +  (l-b3)p3  =  a131  +  a232  »  (3.13) 

b3(l-b2)p2  +  b2(l-b3)p3  =  (a]  +  a2)6132  .  (3.14) 

A  difficulty  with  this  apparently  straightforward  solution  is  that  the  inversion 
of  (3.9)  or  (3.10)  could  lead  to  a  function  which  is  not  a  probability  density, 
or  it  could  yield  a  probability  density  but  not  one  which  is  a  probabilistic 
mixture  of  three  exponentials.   In  fact,  neither  of  these  possibilities  is 
the  case,  as  will  be  shown  by  establishing  that  p2  and  p3  are  positive 
and  subject  to  the  condition  p2  +  P3  <  1  ,  and  hence  can  represent 
probabilities. 


Explicit  expressions  for  p2  and  p3  can  be  obtained  from  (3.13)  and 
(3.14)  and  are  given  at  (3.3)  and  (3.4).  From  now  on  it  will  be  assumed, 
in  accordance  with  the  theorem,  that  b2  is  the  larger  of  b?     and  b.  , 
these  being  obtained  by  solving  the  quadratic  pair  (3.11)  and  (3.12).  To 
establish  that  p2  +  p3  <  1  ,  we  have,  by  adding  (3.3)  and  (3.4), 

( 06-j  3-|  +  a232)  ~  (a]  +  ap)ft-|$2 
p2  +  P3  =  (l-b2)(1-b3) '  <3'15) 

Multiplying  out  (l-b2)(l-b3)  in  the  denominator  and  using  (3.11)  and  (3.12) 
gives,  after  some  rearrangement, 

(1-6-,  )(l-32) 
P2  +  P3  =  ""  _  (1-e-,)(l-32)  +  a131d-B2)  +  a^d-B^  '  (3J6) 

The  algebraic  expression  here  is  clearly  positive  and  less  than  one,  from 
which  it  follows  that  p-,  +  p2  <  1  . 

The  positivity  of  p2  and  p3  will  now  be  proved  by  showing  that  the 
numerators  and  denominators  of  (3.3)  and  (3.4)  are  positive.  For  the 
denominators,  this  requires  that  0  <  b2 ,  b3  <  1  which  will  be  verified 
by  showing  that  0  <  b2b^  <  1  and  0  <  (l-b2)(l-b3)  <  1  .  The  first  of 
these  latter  two  inequalities  is  obvious  from  (3.12);  for  the  second  consider 
the  expressions 

(l-b2)(l-b3)  =  1  -  (b2  +  b3)  +  b2b3 

=  (a-,6-,  +  l-3-|  )(a2B2  +  1  -  32)  -  (a^-, )  (a^g)        (3.17) 

after  using  (3.11)  and  (3.12),  and  then 
1  -  (l-b2)(l-b3)  =  b2  +  b3  -  b2b3 

=  (l-a1)3-,(l-32)  +  (l-OgjBgd-B^  +  B]  B2  .  (3.18) 
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Both  the  right  hand  sides  of  (3.17)  and  (3.18)  are  obviously  positive.  This 
concludes  the  proof  that  0  <  b-j ,  b2  <  1  and  hence  that  the  denominators  of 
P2  and  p3  are  positive. 

For  the  numerators  of  p2  and  p3  to  be  positive  (3.3)  and  (3.4) 
indicate  that  b  =  (c^  +  a^e^/U^  +  a2^  must  satisfy  the  inequalities 

b3  <  b  <  b2  .  (3.19) 

At  this  last  stage,  explicity  expressions  for  b2  and  b~  must  be  used, 
and  from  (3.11)  and  (3.12)  are  given,  after  writing 

S   =    (l-a-jje-j    +   (l-a2)e2      and     r  =    (l~ci-|-a2)3-|62    , 
by 

b2  =  ^js+  (s2-4r)1/2[  and     b3  =  l)s  -   ( s2-4r)1/2[ .  (3.20) 

Then  (3.19)   is  equivalent  to 

-  (s2-4r)1/2  <    s  -  2b  <  (s2-4r)1/2 
or  s2-4r  >   (s-2b)2 

or  sb  -  b2  -  r  >  0  .  (3.21) 

After  some  algebraic  rearrangement  the  left  hand  side  of  (3.21)  becomes 


a1a23132(31-32)2/(a131+a232)2  (3.22) 


which  is  again  clearly  strictly  positive,  excluding  special  cases,  as 
was  to  be  proved. 
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This  concludes  the  proof  that  p2  and  p3  are  both  positive  subject  to 
p2  +  P3  <  1  ;  hence  I-P2  -  P3  ,  Po  ,  and  p^  can  all  be  regarded  as  prob- 
abilities. Thus  E   has  a  proper  probability  distribution  which  can  be 
generated  as  the  (l-p2  -  P3>  P2»  P3)  mixture  of  three  exponentials  of 
means  1,  b?  and  b~  respectively;  further,  both  b2  and  b^  are  less  than 
unity  and  b2  f   b3  . 

In  the  special  cases  mentioned  earlier,  there  are  valid  and  simpler 
results  for  the  distribution  of  E  .  For  instance,  when  3,  =  32  =  1  , 
E   has  a  simple  exponential  distribution  of  mean  (1-a,  -  a«)  .  When 
3,  =  62  f   1  the  innovation  has  a  mixed  exponential  distribution  of  the 
NEAR(l)  form  given  in  (2.1)  with  a  =  a,  +  a2  .  Delineation  of  all  these 
special  cases  is  needed  for  successful  simulation  of  the  NEAR(2)  process 
for  the  full  parameter  range.  This  will  be  considered  elsewhere. 
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4.  Other  Uses  of  the  NEAR(2)  Exponential  Construction 

The  NEAR(2)  process  was  established  by  showing  that  (3.2)  was  a  valid 
innovation  distribution  for  the  relation  (3.1).  The  distributional  assump- 
tions implied  by  this  result  can  also  be  taken  out  of  the  time  series  context  in 
which  they  were  derived  and  viewed  generally  as  a  way  to  combine  a  pair  of 
(possibly  dependent)  unit  exponential  variables  (L, ,  l_2)  with  an  independent 
triple  of  (possibly  dependent)  unit  exponential  variables  (M, ,  M~,  M,J  so 
as  to  yield  a  further  unit  exponential  variable.  Specifically,  with 
(a-,,  a2,  3i»  3o)  and  (bo»  b3'  p2'  p3^  as  Previously  related  by  (3.3)  -  (3.7), 
the  Theorem  has  established  that 


i2L2 


w.p. 
w.p. 
w.p. 


ul 

1-a 


-  a. 


•     M, 

w.p. 

1-P1 

■  +  ■ 

b2M2 

w.p. 

p2 

b3M3 

w.p. 

p3 

(4.1) 


has  a  unit  exponential  distribution. 

First  of  all,  the  idea  of  "switching"  will  be  illustrated;  in  the  NEAR(2) 
context,  this  suggests  taking  (M, ,  M2>  M3)  as  (X  , ,  X  ?,  xn-3^  and 
(L, ,  L2)  as  (E  ,  E  )  .  Then  (4.1)  gives  the  time  series  model 


X..  = 


Vi 

1b2Xn-2 
b3Xn-3 


w.p.  l-p2  -  p3 
w.p.  p2 
w.p.  p3 


H 


!lEn  w'p-  al 
52En  w.p.  a2 


W.p.  1-a-,  -  a2 


(4.2) 


This  is  a  third-order  autoregression,  actually  a  case  of  the  EAR(3)  model 
cited  in  Lawrance  and  Lewis  (1979);  note,  however,  that  this  third-order 
autoregressive  exponential  process  allows  zero  innovations.  Another,  better 
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behaved  higher-order  exponential  model  -in  fact  a  p-th  order  model-  is 
obtained  by  the  following  application  of  the  result  (4.1)  in  its  original 
form  (3.1).  Let  the  indices  l,2,...,p  be  partitioned  into  two  non-empty 
sets  I,  and  I~  of  size  t-,  and  t2  respectively.  Then  in  the  model 


fs^x 


Xn  =  ^ 


1V1  w-p-    al 


*pVp 


w.p 


W.p.    l-oti-...a 


>   +  En  ,   n  =  0,  +  1,  +  2,...,   (4.3) 


Pi 


let  &\  =   3r  lelj  ;  3^-  =  S2 


,  i€l?;  £   a!  -  a-,  and   Y   cu  =  a9 

c   i€i.  ]   '    ui2  '  l 


Then 


if  a^  +  a^  <  1  ,  0  <  3-j ,  32  <  ^  tne  distribution  of  E   is  given  by  the 
Theorem.  Thus  we  have  a  pth-order  exponential  autoregressive  process  with 
four-parameters.  However,  while  this  may  seem  satisfying  it  is  not  clear 
that  four  parameters  would  be  sufficient  to  characterize  the  sample  path 
behaviour  of  an  exponential  process  with  yery   high  order  dependence. 

Another  use  of  (4.1)  is  to  allow  L-,  and  l_2  to  both  be  X  ,  ,  and 
so  obtain  a  four  parameter  first-order  model  of  the  form 


V 


Vn-1 

w.p. 

al 

r  e 

n 

w.p. 

32Xn-l 

w.p. 

a^ 

+  > 

b0E 
2  n 

w.p. 

0 

w.p. 

l-a-.    ■ 

-   a~ 

b0E 
3  n 

w.p. 

l-p9  -  p 


'2  » 


3  » 


(4.4) 


Four  parameters  seems  excessive  for  a  first-order  autoregressive  process  but 
simulations  show  a  wide  range  of  behaviour  in  sample  paths  with  different 
choices  of  parameters.  Equation  (4.4)  in  turn  suggests  a  first-order  model 
allowing  negative  dependence.  This  is  obtained  by  letting  the  variable  in 
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(4.4)  which  is  multiplied  by  e9  be  the  antithetic  transformation  of  X  .  , 

"Vi  n 

that  is  -  log(l-e    )  .  A  two  parameter  version  of  the  model  could  be 
obtained,  for  example,  taking  a,  =  a«  ,  Bi  =  3?  . 

A  third  type  of  use  of  the  construction  is  to  give  mixed  autoregressive 
moving  average  models;  for  this,  (L,,  l_2)  is  (X  ,  ,  X  2)  as  previously, 
but  (M-|,  M2,  M3)  are  chosen  to  be  (E   En+1 »  En+-| )  for  a  second-order  moving 
average  compoent,  or  as  (E  ,  E   ,  En+2^  ^or  a  third-order  moving  average 
component. 

Out  of  the  time  series  context,  the  construction  suggests  ways  to  obtain 
multivariate  exponential  distributions,  rather  as  in  Lawrance  and  Lewis 
(1983b). 

Further  possibilities  are  numerous,  but  it  is  not  the  intention  here  to 
exhaustively  list  them,  or  to  derive  the  details  of  those  cited  at  this  time. 
Analysis  in  the  following  sections  will  deal  with  the  basic  NEAR(2)  model. 
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5.  Autocorrelation  Structure  of  the  NEAR(2)  process 

In  this  section  we  show  that  the  autocorrelations  p(z)   =  corr(X  ,X   ), 

v  n'  n-£  ' 

i  =   0,  +  1,  +  2,...  of  the  NEAR(2)  process  satisfy  Yule-Walker  type  dif- 
ference equations;  thus  the  second  moment  dependency  aspects  of  the  process 
are  indistinguishable  from  those  of  a  standard  autoregressive  model,  AR(2). 
To  show  this  it  is  convenient  to  write  the  equation  (3.1)  as  an  additive 
autoregressive  combination  of  X  , ,  X  «  anc'  En  •  Thus  we  have  a  random 
coefficient,  linear  additive  autoregressive  process 


where 


Xn  '   BlKnXn-l   +  B2KnXn.2  ♦  L|1En 

1   w.p.  l-p2  -  p3  , 

< b2  w.p.  p2       , 

b3  w.p.  p3       ; 


n  =  0,  +  1,  +  2,..., 


n  =  0,  +  1,  +  2,..., 


(5.1) 


(5.2) 


{Kn'  V 


(1,0)     W.p.   a1  , 
•((0,1)    W.p.   a2  » 
(0,0)    W.p.   1-a-j  -  a2  ; 


n  =  0,  +  1 ,  +  2, 


(5.3) 


and  the  i.i.d.  sequences  {L  }  and  {K  ,  K  }  are  mutually  independent  and 

independent  of  the  independent  exponential  sequence  {E  }  .  The  E  '  s  are 

assumed  to  have  mean  and  variance  1  ,  as  do  the  X  's  by  construction, 
i  i 

Now  E(K  )  =  a,  and  E(K  )  =  a2  ,  so  that  E(LR)  =  1-pt^^  -  a262  . 

Then  multiplying  Xn  in  (5.1)  by  Xn_^  we  have,  for  I  >   1  , 


E<Vn-.e>  ■  W<VlViJ  +  a2B2E<Xn-2W  +  ^nl^n^'W 
=  »l6lE(Xn-lW  +  a2B2E(Xn-2Xn-£>  +  ]  "  alBl  "  a2B2  • 
so  that 

E<XnXn-£>  "  ]  "  «lBl{E<Xn-lXn-£>  "  1}  +  "2B21E(Xn-2Xn-Z>   "  1} 
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and 


p(-1)   =  alP(-U-l))  +  a2p(-tf-2))  , 


where  a-,  =  a^  and  a2  =  a^  .  Using  p(-i)   =  p(i) ,  i  >  ^    ,  we  have 
finally 

pU)  =  a]PU-l)  +  a2p(£-2)  ,     I  >  ]    ,  (5.4) 

which  are  the  Yule-Walker  equations  for  the  AR(2)  process.  The  conditions 
for  a  solution  to  exist  (Box  and  Jenkins,  1970,  p. 58),  a,  +  a?  <  1 
a1  -  a2  >  -  1  ,  a2  >  -  1  are  clearly  satisfied  when  the  conditions  on 
c*-j»  «2>  P|.  32  given  in  the  Theorem  of  Section  3  hold.  We  then  have 

p(1)  =  a^O-a^  and  P(2)  =  a^O)  +  a2  .  (5.5) 

Note,  however,  the  restriction  to  positive  correlations  since  a,  and  a? 
are  positive.  The  possible  region  of  (p(l),p(2))  values  is  bounded  below  by 
p(2)  =  p  (1)  and  otherwise  bounded  by  p(l)  >  0  and  p(2)  <_  1  .  Broadening 
of  the  model  to  negative  dependency  may  be  achieved  using  antithetic  ideas, 
or  the  bivariate  scheme  given  in  Gaver  and  Lewis  (1980),  but  is  not  pursued 
here. 

Note  too  that  the  parameters  in  (5.4)  enter  only  as  products  a,  =  a-,6-1 
and  a2  =  a23?  .  Thus  for  small  enough  a,  and  a2  ,  values  of  B-,  and  b^ 
greater  than  one  could  be  allowed,  and  (5.4)  would  still  have  a  stable 
solution.  However  the  sequence  E   in  the  defining  equation  (3.2)  may  not 
exist;  it  is  not  known  if  8,  <_  1  and  32  <_  1  is  a  necessary  condition  for 
this  existence. 

Solutions  to  (5.4)  are  given  in  detail  in  Box  and  Jenkins  (1970, 
pp.  58-60). 

17 


Specifying  allowable  values  of  p(l)  and  p(2),  as  may  be  done  by  an 
initial  second  order  analysis  of  data,  leaves  two  parameters  to  be  specified 
in  the  model,  say  a-,  and  a2  which  could  produce  very   different  sample 
path  behaviour  in  the  time  series.  It  is  important  to  notice  that  this 
specification  of  p(l)  and  p(2)  further  constrains  the  range  of  possible 
a,  and  a2  values.  Recalling  that  p(l)  and  p(2)  fix  a-,  =  di3i  and 
a?  =  aoBo'  as  well  as  that  a,  +  a2  £  1  ,  it  is  easily  shown  that  we  must 
have 


a-,  <_  a-i  and  a~  <_  a~ 


which  implies  that  a-,  +  a?  <_  a-,  +  a2  <_  1  .  Thus  a1  and  a?  are  forced 
to  lie  in  a  triangular  subregion  of  the  triangular  (a, ,a2)  region  which  is 
bounded  below  by  a?  ,  bounded  on  the  left  by  a-,  ,  and  bounded  above  by 
the  1  ine  a-,   +  ou  =  1  • 
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6.  An  analysis  of  wind  velocity  data 
6.1.  Discussion  of  the  data 

Lewis  and  Hugus  (1983)  have  given  an  analysis  of  a  set  of  3-hourly 
wind  velocity  readings  taken  by  ship  PAPA  in  the  Gulf  of  Alaska  over  a  15 
year  period.  After  suitable  detrending  to  remove  1  year,  6  month,  12  hour 
and  6  hour  cyclic  trend  components  a  first-order  autoregressive  Gamma  model 
(Lewis,  1981)  was  fitted  to  the  data,  the  use  of  the  model  being  suggested 
by  the  shape  of  the  (marginal)  histogram  of  the  data  (Figure  6.1)  and  the  shape 
of  the  normalized  log  periodogram  of  the  data  (Figure  6.2).  Note  that  the 
sample  size  is  N  =  43,800;  also  there  is  a  residual  6-hour  effect  (P  =  21900) 
because  this  cycle  varies  in  intensity  over  the  15  years.  In  what  follows 
this  will  be  ignored  and  the  data  will  be  treated  as  stationary. 

It  is  not  the  object  here  to  discuss  the  above  analysis  in  detail  but 
to  discuss  a  different  analysis  of  the  data  using  an  assumption  of  a  Weibull 
marginal  distribution  and  a  transformation  to  exponential  variables.  This 
is  suggested  firstly  by  the  fact  that  a  Weibull  distribution  is  commonly  used 
for  this  type  of  data  by  meteorologists  and  secondly  by  the  fact  that  Weibull 

and  Gamma  distributions  fit  the  data  equally  well  (Lewis  and  Hugus,  1983). 

'     2  ?f)5 
The  histogram  of  the  transformed  data,  X  =  X  '   ,  is  shown  in  Figure  6,3, 

where  the  power  transformation  to  exponential ity  has  been  determined  by  fitting 

the  empirical  coefficient  of  variation,  0.479,  to  the  theoretical  Weibull 

coefficient  of  variation,  C(X)  =  {r(2/c  +  l)/[r(l/c  +  I)]2}  -  1  to  give 

c  =  2.205.  This  transformation  does  affect  the  correlation  structure  of  the 

data,  as  shown  in  Table  6.1. 

Column  2  in  Table  6.1  gives  the  estimated  auto  correlations,  pU)  , 

of  the  detrended  data;  rough  95%  confidence  bands  for  these  estimates  are 

given  by  adding  and  subtracting  2/(N)1/2  =  2/(43, 800)1/2  ~  .01  .  The  first 
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column  gives  the  fitted  auto  correlations  for  a  model  with  AR(l)-type  auto 

correlations,  just  pU)  =  (p(l))£  =  (.8214)*  ,  for  lags  i  =   1,.  ..,5. 

The  maximum  difference  (3rd  column)  is  0.0109  and  the  fit  is  clearly  good. 

Column  6  in  Table  6.1  gives  the  estimated  auto  corelations,  p'U)  , 
for  the  transformed  data;  the  transformation  changes  the  serial  correlations 
slightly.  However,  columns  5  and  7,  which  give  the  fitted  AR(1)  and  AR(2) 
correlation  models  respectively,  show  that  a  model  with  AR(1)  correlation 
structure  is  not  adequate,  but  that  a  model  with  AR(2)  correlation  fits  well. 
This  fit  is  borne  out  by  a  periodogram  analysis,  which  is  not  given  here. 

Thus,  a  NEAR(2)  model  is  a  candidate  for  representing  the  transformed 
data,  and  if  p(l)  and  p(2)  are  fixed  at  the  estimated  values  of  0.7985 
and  0.6574,  then  the  corresponding  a-,  =  (a-,3-,)  and  a2  =  (a232)  "from  (5.5) 
are,  respectively,  0.75488  and  0.05463.  There  are  still  two  degrees  of  free- 
dom left  in  fitting  the  model,  represented  by  choice  of  parameters  a,  and 
a~  greater  than  an  equal  to  0.75488  and  0.05463  respectively,  with 

°h  +  a2  —  ^  • 

Figure  6.4  shows  the  logarithm  of  the  periodogram  of  the  usual  AR(2) 

model  linear  residuals,  R  =  X'  -  anX'  -,  -  a0X'   0,  of  the  transformed  data. 

n    n    I  n- 1    I  x\-l 

This  is  flat  (ignoring  the  slight  effect  at  period  6  hours).  At  this  point 
it  might  be  thought  that  the  usual  second-order  autoregressive  model  is 
adequate.  We  shall  however  now  develop  an  extended  residual  analysis  for 
higher  order  dependence  which  justifies  fitting  the  NEAR(2)  model  to  the 
transformed  data. 
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7.  A  Residual  Analysis  for  the  NEAR(2)  Model 
7.1.  General  results 

It  has  already  been  remarked  that  the  autocorrelations  p(jt)  are 
insufficient  to  describe  the  dependency  structure  of  NEAR(2)  models.  A 
natural  next  step  might  be  to  examine  higher  joint  moments  and  consider  their 
associated  spectra  (see  e.g.,  Priestley,  1982).  The  functions  so  obtained,  e.g. 
the  bispectrum,  are  often  found  to  be  difficult  to  calculate  and  hard  to 
interpret.  Rather  than  follow  this  course,  it  is  proposed  to  adapt  some 
ideas  from  a  residual  analysis  for  autoregressive  models  which  has  recently 
been  suggested  by  the  authors  (Lewis  and  Lawrance,  1983).  The  thrust  of 
this  analysis  is  that  the  standard  process  of  fitting  and  validating  a  linear 
autoregressive  model  is  carried  out  beyond  the  customary  final  stage  at  which 
uncorrelated  residuals  are  obtained  (as  in  the  previous  section).  The  usual 
presumption  is  that  the  residuals  are  not  only  uncorrelated  but  also  inde- 
pendent. This  need  not  be  the  case, as  will  be  exemplified  for  the  wind  ve- 
locity data.  Moreover,  dependent  but  uncorrelated  residuals  are  obtained 
(Lewis  and  Lawrance,  1983)  even  for  the  NEAR(2)  process.  Thus  the  residuals 
should  be  subjected  to  further  analysis  in  respect  of  their  remaining  de- 
pendency. Any  found  is  then  evidence  that  a  standard  linear,  constant  co- 
efficient second-order  autoregressive  model  is  deficient.  With  normally 
distributed  time  series  data  this  might  suggest  that  non-linear  modelling 
should  be  explored.  With  data  marginally  distributed  in  some  other  identi- 
fiable manner,  then  the  exploration  of  a  particular  type  of  model  with 
specified  marginal  distribution  and  autocorrelation  function  is  suggested. 
This  latter  course  is  envisaged  here. 
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Higher  order  dependency  properties  of  the  uncorrected  residuals  are 
obtained  for  the  proposed  model  and  compared  with  their  data  counterparts; 
this  stage  can  be  informative  from  both  exploratory  and  estimation  consid- 
erations, and  can  be  thought  of  as  part  of  the  model -refinement  process 
common  to  much  statistical  methodology. 

It  might  be  thought  that  the  specific  class  of  NEAR(2)  models  could  be 
incorporated  in  a  residual  analysis  in  the  standard  manner.  However,  a 
moments  reflection  indicates  that  even  after  estimating  parameters  of  (2.2), 
it  will  not  be  possible, because  of  the  mixture  involved,  to  write  down  an 
expression  or  recursion  for  the  residual  innovations,  namely  E     .  However, 
the  corresponding  autoregressive  (or  linearized)  residual  is  available  and 
given,  as  in  the  previous  section,  by 

Rn  =  Xn  "  alVl  "  a2Xn-2  •  (7'1> 

We  now  show  that  these  are  uncorrelated  for  the  NEAR(2)  process. 


7.2.  The  residual  theorem. 

The  autocovariances  of  the  residuals  (7.1)  are 

Cov(Rn,  Rn+£)  ■  Cov(XR,  Rn+£)  -  a1Cov(Xn_1,  RR+£)  -  a2Cov(Xn_2,  Rn+£) 

=  Cov(V  W  "  alCov(V  Wl}  -  a2Cov(V  W2}  ' 

since  the  {X  }  process  and  consequently  the  {Rn}  process  is  stationary. 
The  covariances  on  the  right  hand  side  are  all  of  the  same  type  and  given  by 


(7.2) 


Cov(Xn,  Rn+£)  =  Cov{Xn,  (Xn+£  -  a^^  -  a2Xn+£_2)} 

(7.3) 

=  {Var(X)}(p(£)  -  a]P(£-l)  -  a2p(£-2)),  z  =   1,2,... 


^7 


By  the  Yule-Walker  equations  (5.4)  the  expression  in  brackets  is  zero,  and 
hence  also 

Corr(Rn,  Rn+£)  =  0  ,    £=+l,+2,...  (7.4) 

as  was  to  be  proved.  That  these  residuals  are  uncorrelated  is  an  immediate 
consequence  of  the  autocorrelations  following  Yule-Walker  equations;  this 
emphasizes  that  this  type  of  residuals  will  be  uncorrelated  for  any  model 
whose  autocorrelations  satisfy  Yule-Walker  equations. 

The  analysis  of  the  uncorrelated  residuals  R  ,n  =  3,4,...  should  begin 
with  scatter  plots  of  the  low  lag  adjacent  values;  any  patterns  or  concen- 
trations will  be  evidence  of  dependency  in  the  residuals.  Many  possibilities 
present  themselves  but  only  one  is  pursued  in  the  following  Section  8,  and 
then  applied  in  Section  9  to  a  continued  analysis  of  the  wind  velocity  data. 


28 


2 

8.  Cross-correlation  analysis  of  {R  }  and  {R  } 

After  the  satisfactory  fit  to  data  of  an  ordinary  linear  model,  the 
residual,  Rn>  should  not  only  be  uncorrelated  but  also  independent;  the 
latter  is  customarily  investigated  by  seeking  a  flat  spectrum,  while  for 
the  independence,  a  flat  spectrum  of  the  squared  residuals  can  be  sought. 
As  a  method  for  probing  model  validity,  the  examination  of  squared  residuals 
has  been  employed  by  McLeod  and  Li  (1983),  following  Cranger  and  Andersen 
(1978);  these  authors  suggested  bilinear  modelling  for  dependent  but  un- 
correlated residuals  of  ARMA  models.  It  is  suggested  here  that  autocorrela- 
tions of  squared  residuals  and  cross  correlations  of  residuals  of  squared 
residuals  be  used  in  the  analysis  of  higher  order  dependence  of  the  detrended 
transformed  wind  speed  data  and  its  modelling  by  the  NEAR(2)  process.  The 
residuals  for  this  data  have  already  been  shown  to  have  a  flat  spectrum 
(Figure  6.4)  while  the  curved  plot  in  Figure  8.1  for  the  cumulative  periodogram 
shows  that  the  spectrum  of  the  squared  residuals  is  far  from  flat. 

Theoretical  investigation  of  the  squared  residuals  of  the  NEAR(2)  model 
is  pursued  here.  Whilst  the  autocorrelations  of  the  squared  residuals  have 
just  been  mentioned,  for  the  NEAR(2)  model  this  involves  computation  of 
36  terms,  mostly  distinct  types  of  ^-  order  moments.  A  simpler  suggestion 
which  involves  only  third-order  moments,  and  is  thus  the  next  step  up  after 
autocorrelations,  is  to  use  the  cross  correlations  of  the  Rp  and  Rn  sequences; 
apart  from  lag  0,  zero  values  will  be  found  for  linear  models.  This  sug- 
gestion is  more  tractable  than  the  squared  residual  analysis  and  will  now 
be  described.  Sampling  properties  of  third-order  moments  are  also  likely 
to  be  less  extreme  than  those  of  fourth-order  moments. 

The  starting  point  of  the  calculation  is  to  note  that  from  the  definition 
of  Rn  at  (7.1)  that 
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E<RnW  ■  E<RnV*)  "  ¥(<Vm)  "  *2E<RnXn-:-l>  ■        (3-5) 
whence  there  is  the  structural  form, 

Cov(R^,  Rn.()  -  Cov(R2,  Xn.p)  -  alCov(R2,  X^,,)  -  .^(rZ,  x^^)  .   (8.6) 

Calculation  of  the  covariance  terms  in  (8.6)  requires  the  expanding  out  of 

2 
Rn,  taking  expectations  and  expression  in  covariance  terms.  Thus 

RnVe   =  Vn-y,  +  alVlVi  +  a2Xn-2Xn-Ji 

-  2Wn-lVi  -  2a2XnXn-2Xn-e 

+  2ala2Xn-lXn-2Xn-,-  (8.7) 


The  conversion    to  covariances  yields 


Cov(R2,  Xn_£)  =  JjU)  +  afj1(A-l)  +  a^U-2) 


where 


-2a.,J2U)  -  2a2J3(0  +  2a-|a2J2^.-l )    ,  (8.3) 


^(0   =  Cov(X2,  Xn_A)    ;     J2U)   =  Cov(XnXn_r  Xp ,,)    . 

J3(e)   =  Cov(XnXn  2,   Xp_£)    .  (3.9) 


We  thus  see  the  types  of  third-order  joint  moments  which  are  involved  in 

2 

the  Cov(R  ,  R   )  calculation:  each  of  these  has  to  be  obtained  using  two 

difference  equations,  one  for  l  >   0  and  one  for  ?  <  0  . 

Taking  J-.U)  for  illustration,  square  each  side  of  the  NEAR(2)  defin- 
ing equation  (2.2)  and  multiply  by  X  .  After  converting  to  the  required 
covariances,  the  recursion  is  found  to  be 

J1  («.)  =  a13]J1(£-l)  +  a^jy*1^)  +  2(l-a1  -  a2)pU)  .  (8.10) 
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For  negative  £  there  is  the  corresponding  equation 


J^U)  =  a1J1  (£+1)  +  a^U+2)  , 


8.11) 


obtained  by  multiplying  each  side  of  (2.2)  by  X    and  converting  to 
covariances.  These  equations  are  given  for  illustration:  there  are  similar 

equations  for  J~{  a)     and  J3U)  ,  and  various  special  cases.  Also  required 

2 
is  the  variance  of  R   which  involves  many  terms. 

The  complete  algorithm  for  computing  the  cross-covariances  is  as  follows. 


Algorithm  for  computing  cross-covariances  for  £  =  0,  +  1,  ...  +  L  . 
0.  Input  a, ,  a,,,  g, ,  32  • 

Set  p(0)  =  1;  p(1)  =  a^/Cl-ag);  p(2)  =  a-jp-j  +  a2  and,  for  £  =  3,...,K  , 
p( jl)  =  a,p(£-l)  +  a2p(£-2).  Note  that  p(-£)  =  pU)  but  these  are 
not  needed  in  computations. 


1.  Set  jyo)  =  4;  J(l)  =  [Aa^^+a^)  +   2(l-a1-a2)]/(l-a2e2)  . 


2a.     J,U)  =  a-,  J-,  ( £+1 )  +  a2J1  ( ji+2 
2b.     J-j  ( £,)  =  a-|  e-|  J-j  ( 52.-1 )  +  a2e2J1 


3.     Set  J2(0)   =  J^l)  -  p(l)  +  1 


4a.     J2U)  =  a-j  J2(£+l )  +  a2J2(£+2 

4b.     J2(£)  =  a-|J,(fc-l)  +  a2J2(£-l 

-[1+p(D]  +   [1  +  pU 

5.     Set  J3(0)  =  ^(2)  -  p(2)  +  1 

6a.     J3U)  =  a1J3(£+l)  +  a2J3(£+2 

6b.     J3(£)  =  a-jJ2(£-l)  +  a2J1(£-2 

+    [l+P(£-2)](l-ara2 


for         £  =  -1,  -2,...,  -   (L+4)    . 
£-2)  +  2(l-a1-a2)P(£)         for         £=2,3,...L 


J2(D  =  jy-D  -  p(D  +  1 


for         £  =  -1,-2,...,  -(L-3)    . 
+  2a1   +  (l+p(l))a2 
l)](l-a.,-a2)         for         £  =  2,3,...,  L   . 

J3(l)   =  J2(2)+p(D-p(2);   J3(2)   =  J1(-2)-p(2)  + 

for         £  =  -1,  -2,...,   -   (L+2)    . 
+    [l+p(l)]a1   +  2a2 
-   (l+p(2))         for         £  =  3,4,5,...,   L   . 
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7.  JU)  =  J-[{a)  +   a^U-1)  +  apii-2)   -  2a]J2(£) 

-  2a2J3U)  +  2a1a2J2(£-l)    for    a   =  -  (L+2),...,  (L) 

2 
8.  Cross-covariance  (Rn,  R   )  =  J(n)  -  a^U+1)  -  a2J(£+2) 

for    £  =  -  L,. . . ,  L  . 
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9.  Further  Analysis  and  Modelling  of  the  Wind  Velocity  Data 

Dependence  in  the  uncorrelated  linear  residuals  for  the  NEAR(2)  model 

has  already  been  demonstrated  (Figure  8.1).  Further  evidence  of  this  is 

2 
provided  by  the  cross-correlations  of  R   and  R   given  in  Figure  9.1. 

The  corresponding  theoretical  cross-correlations  for  the  NEAR(2)  model  will 

next  be  presented,  having  been  computed  using  the  algorithm  given  in 

Section  8.  At  this  point,  it  will  be  recalled,  the  NEAR(2)  model  has  not 

been  fitted  in  terms  of  all  4  parameter;  the  residuals  involve  the  model 

parameter  only  through  a-.  =  a-,3-1  and  a2  =  a232  and  a-., 6-j  »««,&«  have 

not  been  separately  estimated.  In  the  present  rather  exploratory  analysis 

the  estimation  problem  will  be  circumvented;  the  cross-correlations  of 

(R  ,R  )  will  be  given  for  four  representative  sets  of  parameter  values  in 

the  reduced  allowable  region,  as  constrained  by  (a-,  >_  a-, ,  a2  >  a2,  a,  +  a«  <.  1) 

according  to  (5.6).  For  the  detrended  transformed  wind  speed  data, 

a-,  =  0.75488  at  a2  =  0.05463,  and  the  four  chosen  sets  are  given  by 


(A)  a,  =  0.760,  a2  =  0.06 

(B)  a,  =  0.925,  a2  =  0.06 

(C)  a-,  =  0.760,  a2  =  0.20 

(D)  a,  =  0.850,  a2  =  0.10 


3,  =  0.99326,  32  =  0.91050 

3n  =  0.81608,  32  =  0.91050 

31  =  0.99326,  32  =  0.27315 

31  =  0.88809,  32  =  0.54636. 


2 
In  Figure  9.2  the  cross-correlations  of  (R  ,R  )  for  each  of  these 

cases  are  presented;  there  is  considerable  differentiation  amongst  the 

figures,  although  the  zero  values  for  all  but  one  or  two  negative  lags  are 

evident  in  all  four  cases.  This  effect  will  be  due  to  the  directional 

nature  of  the  process.  For  matching  with  the  wind  speed  data,  Case  D 

appear  the  most  promising,  although  the  suggestion  of  a  negative  correlation 

at  lag  minus  one  lends  a  little  support  to  Case  A  .  Cases  B  and  C  are 
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definitely  unsympathetic  to  the  data.  Hence,  a  choice  of  parameters  inter- 
mediate between  Cases  A  and  D  is  suggested  by  this  exploratory  analysis. 
A  fuller  analysis  would  require  estimates  of  all  four  parameters  and  compari 
son  of  the  resulting  residual-squared  residual  cross-correlations  with  the 
estimated  cross-correlations  obtained  from  the  data,  as  shown  in  Figure  9.1. 
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10.  Conclusions  and  Further  Analysis 

The  very   broad  four  parameter  NEAR(2)  time  series  model  having  expo- 
nential marginals  and  the  correlation  structure  of  an  AR(2)  model  has  been 
established;  further  developments  will  be  reported  elsewhere.  A  preliminary 
fit  of  the  NEAR(2)  model  has  been  made  to  a  very  long  series  of  wind  speed 
data,  the  data  having  been  detrended  and  transformed  so  as  to  have  exponen- 
tially distributed  marginals;  utility  of  the  suggested  residual  analysis  in 
probing  higher  order  dependence  has  been  demonstrated.  This  residual  analysis 
has  been  based  on  the  cross-correlations  between  the  residuals  and  the 
squared-residuals. 

An  extension  to  this  analysis  using  reversed  residuals  is  possible; 
more  of  the  higher  order  dependency  of  the  NEAR(2)  model  would  be  revealed 
and  this  would  enable  further  aspects  of  its  suitability  of  the  model  for 
the  wind  speed  data  to  be  ascertained. 
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